# Установим рабочую директорию
setwd("~/Documents/GitHub/RProjects_EF/Networks")
# Загрузим данные для исследования: набор ребер и вершин
edges <- read.csv("edges_1.csv", header=T)
colnames(edges) <- c('from', 'to', 'type1', 'type2', 'weight')
nodes <- unique(edges[c('from', 'type1')])
more_nodes <- unique(edges[c('to', 'type2')])
colnames(nodes) <- c('organisation', 'type')
colnames(more_nodes) <- c('organisation', 'type')
nodes <- rbind(nodes, more_nodes)
nodes <- unique(nodes[c('organisation', 'type')])
rownames(nodes) <- NULL
edges <- edges[c('from', 'to', 'weight')]
# Создадим также набор обратных величин для весов ребер (это пригодится нам в дальнейшем). Исходные веса нам уже не понадобятся
edges$weight <- 1 / edges$weight
net <- graph_from_data_frame(d=edges, vertices=nodes, directed=F)
V(net)$color <- ifelse( V(net)$type=="Business Enterprise", "red",
ifelse( V(net)$type=="Private not for profit", "purple",
ifelse( V(net)$type=="Government", "green",
ifelse( V(net)$type=="Higher Education", "blue", "black")
)
)
)
head(V(net)$color, 10)
[1] "red" "red" "red" "red" "red" "red" "purple" "green" "purple" "red"
Поясним, почему в дальнейшем будем работать именно с обратными весами: Кратчайший путь и кратчайшее расстояние в данной задаче будут представлять собой не то же самое, что и в задаче с издержками перевозок или длинами расстояний между вершинами в качестве весов ребер. Отрицательные веса мы не используем, потому что алгоритм Дейкстры работает с неотрицательными весами ребер, поэтому далее мы оперируем с обратными весами, т.к., на наш взгляд, этот выбор является наиболее оптимальным для преодоления возникших трудностей.
Исследовательский вопрос: Какие есть способы улучшения взаимодействия учреждений высшего образования Новой Зеландии (если это вообще требуется)?
Наше исследование показывает, что возможно повысить эффективность небольших университетов посредством увеличения сотрудничества с крупными, которые в свою очередь имеют тесное сотредничество с бизнесом.
Детали см. в Задании 6.
Датасет представляет из себя множество вершин с качественными характеристиками и множество ребер с числовыми характеристиками (весами).
Вершинами являются организации Новой Зеландии, которые занимаются написанием научных статей в качестве одного из видов деятельности. Качественной характеристикой вершин является вид организации (государственная организации, высшее учебное заведение, коммерческое предприятие и частное некоммерческое предприятие). Всего в датасете 1511 вершин.
Между i-ой и j-ой вершинами (организациями)
существует ребро, если за период 2010-2015 гг. на Scopus была
опубликована хотя бы одна статья, хотя бы один автор которой числится в
организации i, а другой автор – в организации
j. Всего в датасете 4273 ребра.
Вес ребра равен количеству совместных опубликованных научных работ.
Если вес ребра между i-ой и j-ой вершиной
равен 3, то это означает, что за указанный период на Scopus было
опубликовано 3 научные работы, у которых хотя бы один автор числится в
организации i, а другой — в организации j.
При использовании нескольких пакетов для работы с графами не
оказалось возможным интерпретируемо визуализировать граф полностью.
Поэтому для отображения структуры нашей сети найдем 100 самых
центральных по степени вершин и представим связи между ними при помощи
пакета ggraph.
V(net)$degrees <- degree(net)
nodes$degrees <- V(net)$degrees
for_vis_n <- nodes[order(-nodes$degrees),][1:100, 1:2]
for_vis_e <- edges[(edges$from %in% for_vis_n$organisation & edges$to %in% for_vis_n$organisation),]
for_vis_net <- graph_from_data_frame(d=for_vis_e, vertices=for_vis_n, directed=F)
V(for_vis_net)$degrees <- degree(for_vis_net)
for_vis_net %>%
ggraph(layout='graphopt') +
geom_node_point(aes(color = type, size=round(degrees/20))) +
geom_edge_link(alpha = 0.05) +
theme_void() +
ggtitle("Network of New Zealand organisations")
Как можно увидеть, среди 100 самых центральных вершин с количественной точки зрения преобладают организации государственного сектора, однако самыми крупными здесь являются новозеландские университеты. Первые 5 мест по числу совместных научных работ принадлежат именно им, лидером по данному показателю является Университет Окленда, самый крупный и престижный университет Новой Зеландии.
Рассчитаем базовые параметры нашего графа. Так как весами дуг в данном случае является число совместных работ, то есть “сила связи” двух вершин, будем использовать обратные веса для подсчета центральностей через кратчайшие пути. Это довольно разумно, так как от написания дополнительной общей работы прирост в силе связи двух организаций должен монотонно снижаться, как и разница в пути до соседа.
#Считаем базовые параметры графа
# Density
# The proportion of present edges from all possible ties.
p <- edge_density(net)
# связность сети
no <- igraph::components(net)$no # в пакете statnet такое же название
csize <- igraph::components(net)$csize
csize
[1] 1463 2 2 2 3 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
[21] 2 2 3 2
# Diameter (longest geodesic distance) (если несколько компонент, то диаметр самой большой)
# Note that edge weights are used by default, unless set to NA.
diam <- diameter(net, weights=NA)
diam_w <- diameter(net)
path <- get_diameter(net, weights=NA)
path_w <- get_diameter(net)
path
+ 7/1511 vertices, named, from 83998d9:
[1] Dairy Production Systems Ltd. Straven Road Veterinary Centre AgResearch
[4] Auckland Council Auckland District Health Board New Zealand National Eye Bank
[7] Project Litefoot Trust
path_w
+ 9/1511 vertices, named, from 83998d9:
[1] Dairy Production Systems Ltd. Straven Road Veterinary Centre AgResearch
[4] Massey University University of Auckland Victoria University of Wellington
[7] GNS Science ECoast Limited Dumpark Data Science
# рассчитаем кратчайшие пути
# Average path length
# The mean of the shortest distance between each pair of nodes in the network
# (in both directions for directed graphs).
m_dist <- mean_distance(net, weights=NA)
m_dist_w <- mean_distance(net)
Плотность нашего графа составляет 0.0037. Такая низкая плотность (вероятность образования связи) довольно естественная для вероятности написания научной работы. Наша сеть имеет 24 компонент связности: одна главная, соединяющая практически все вершины, а все остальные имеют по 2-3 вершины. Диаметр (наибольший из кратчаших путей) имеет длину 4.01, или 6, если брать единичные веса для всех дуг.
Изучим соседство самой центральной вершины — Оклендского Университета (с кем чаще всего взаимодействует наиболее престижный ВУЗ страны):
# We can also easily identify the immediate neighbors of a vertex.
# The 'neighbors' function finds all nodes one step out from the focal actor.
# To find the neighbors for multiple nodes, use 'adjacent_vertices()'.
d <- data.frame(neighbors(net, V(net)["University of Auckland"], mode="total")$type)
groups <- as.data.frame(table(d))[, 2]
pie(groups, labels=paste0(as.data.frame(table(d))[, 1], ', ',round(groups/sum(groups)*100, digits=1), '%'),
col=c(5, 2, 3, 4), main='Neighbors of University of Auckland')
# To find node neighborhoods going more than one step out, use function 'ego()'
# with parameter 'order' set to the number of steps out to go from the focal node(s).
neigh2 <- data.frame(ego(net, order = 2, nodes = V(net)["University of Auckland"])[[1]]$type)
groups2 <- as.data.frame(table(neigh2))[, 2]
pie(groups2, labels=paste0(as.data.frame(table(neigh2))[, 1], ', ',round(groups2/sum(groups2)*100, digits=1), '%'),
col=c(5, 2, 3, 4), main='Two-step neighbors of University of Auckland')
Довольно ожидаемо, что чаще всего в написании научных работ такой успешный университет сотрудничает с индустрией. Это подтвеждает наше предположение о том, успех научной работы зависит от работы с бизнесом. В ближайшем окружении данной организации также принято вести исследования с индустрией.
Сравним нашу сеть с теоритическими моделями случайного графа на основе распределения числа их соседей, диаметра, среднего кратчайшего пути и числа компонент связности. Начнем с модели случайного графа Эрдеша-Реньи.
#Генерируем модели случайного графа erdos.renyi
p <- edge_density(net)
random_model <- erdos.renyi.game(n=gorder(net), p.or.m=p, type='gnp')
par(mfrow = c(1, 2))
hist(degree(net), main = "Histogram of node degree", breaks = 1:vcount(net)-1, prob = TRUE, xlim=range(c(0, 20)))
hist(degree(random_model), main = "Histogram of node degree", breaks = 1:vcount(net)-1, prob = TRUE, xlim=range(c(0, 20)))
par(mfrow = c(1, 1))
gl <- vector('list', 100)
m <- rep(NA, 100)
d <- rep(NA, 100)
a <- rep(NA, 100)
for (i in 1:100){
gl[[i]] <- erdos.renyi.game(n=gorder(net), p.or.m=p, type='gnp')
m[i] <- mean_distance(gl[[i]])
d[i] <- diameter(gl[[i]])
a[i] <- igraph::components(gl[[i]])$no
}
mean_distance <- m
hist(mean_distance, xlim=range(c(min(m, mean_distance(net)), max(m, mean_distance(net)))))
abline(v=mean_distance(net), col='red', lty=3, lwd=2)
diameter <- d
hist(diameter, xlim=range(c(min(d, diameter(net, weights=NA)), max(d, diameter(net, weights=NA)))))
abline(v=diameter(net, weights=NA), col='red', lty=3, lwd=2)
number_of_components <- a
hist(number_of_components, xlim=range(c(0, 25)))
abline(v=igraph::components(net)$no, col='red', lty=3, lwd=2)
Как можно заметить, все рассматриваемые показатели нашей сети значимо отличаются от модельных, значит, структуру нашей сети нельзя описать данным образом.
Обратимся к модели “Малого мира” с регулярной решеткой. Однако в этом случае сперва нужно подобрать гиперпараметр — долю переприсоединенных дуг.
#Генерируем модели случайного графа с регулярной решеткой
n_of_nei <- mean(degree(net))
diam_net <- diameter(net, weights=NA)
#Ищем оптимальное значение p
mean_dist_net <- mean_distance(net)
model <- vector('list', 100)
m <- rep(NA, 100)
d <- rep(NA, 100)
res <- rep(NA, 100)
for (i in 1:100){
model[[i]] <- sample_smallworld(dim = 1, size = gorder(net), nei = n_of_nei, p = i/100)
m[i] <- mean_distance(model[[i]])
d[i] <- diameter(model[[i]])
res[i] <- ((m[i]-mean_dist_net)^2+(d[i]-diam_net)^2)^(1/2)
}
# p = 0.99
random_model <- sample_smallworld(dim = 1, size = gorder(net), nei = n_of_nei, p = 0.99)
par(mfrow = c(1, 2))
hist(degree(net), main = "Histogram of node degree", breaks = 1:vcount(net)-1, prob = TRUE, xlim=range(c(0, 30)))
hist(degree(random_model), main = "Histogram of node degree", breaks = 1:vcount(net)-1, prob = TRUE, xlim=range(c(0, 30)))
par(mfrow = c(1, 1))
gl <- vector('list', 100)
m <- rep(NA, 100)
d <- rep(NA, 100)
a <- rep(NA, 100)
for (i in 1:100){
gl[[i]] <- sample_smallworld(dim = 1, size = gorder(net), nei = n_of_nei, p = 0.99)
m[i] <- mean_distance(gl[[i]])
d[i] <- diameter(gl[[i]])
a[i] <- igraph::components(gl[[i]])$no
}
mean_distance <- m
hist(mean_distance, xlim=range(c(min(m, mean_distance(net)), max(m, mean_distance(net)))))
abline(v=mean_distance(net), col='red', lty=3, lwd=2)
diameter <- d
hist(diameter, xlim=range(c(min(d, diameter(net, weights=NA)), max(d, diameter(net, weights=NA)))))
abline(v=diameter(net, weights=NA), col='red', lty=3, lwd=2)
number_of_components <- a
hist(number_of_components, xlim=range(c(0, 25)))
abline(v=igraph::components(net)$no, col='red', lty=3, lwd=2)
Аналогично предыдущему случаю, данная модель плохо описывает структуру нашего графа. Значимо совпадает лишь диаметр, так как по нему в том числе происходил подбор гиперпараметра.
Рассмотрим также схожесть с моделью Барабаши-Альберта, включающей предпочтительное присоединение к наиболее центральным вершинам.
# Генерируем модель Барабаши-Альберта
#Ищем оптимальное значение power
model <- vector('list', 100)
res <- rep(NA, 40)
m <- rep(NA, 100)
d <- rep(NA, 100)
for (i in c(1:40)){
model[[i]] <- sample_pa(n = gorder(net), power = i/10, directed = FALSE)
m[i] <- mean_distance(model[[i]])
d[i] <- diameter(model[[i]])
res[i] <- ((m[i]-mean_dist_net)^2+(d[i]-diam_net)^2)^(1/2)
}
#power=2.15
random_model <- sample_pa(n = gorder(net), power = 2.15, directed = FALSE)
par(mfrow = c(1, 2))
hist(degree(net), main = "Histogram of node degree", breaks = 1:vcount(net)-1, prob = TRUE, xlim=range(c(0, 20)))
hist(degree(random_model), main = "Histogram of node degree", breaks = 1:vcount(net)-1, prob = TRUE, xlim=range(c(0, 20)))
par(mfrow = c(1, 1))
gl <- vector('list', 100)
m <- rep(NA, 100)
d <- rep(NA, 100)
a <- rep(NA, 100)
for (i in 1:100){
gl[[i]] <- sample_pa(n = gorder(net), power = 2.15, directed = FALSE)
m[i] <- mean_distance(gl[[i]])
d[i] <- diameter(gl[[i]])
a[i] <- igraph::components(gl[[i]])$no
}
mean_distance <- m
hist(mean_distance, xlim=range(c(min(m, mean_distance(net)), max(m, mean_distance(net)))))
abline(v=mean_distance(net), col='red', lty=3, lwd=2)
diameter <- d
hist(diameter, xlim=range(c(min(d, diameter(net, weights=NA)), max(d, diameter(net, weights=NA)))))
abline(v=diameter(net, weights=NA), col='red', lty=3, lwd=2)
number_of_components <- a
hist(number_of_components, xlim=range(c(0, 25)))
abline(v=igraph::components(net)$no, col='red', lty=3, lwd=2)
Как мы получили, данная модель тоже не подходит для описания структуры новозеландской сети: все статистически значимо отличается, кроме параметра, по которому проходил подбор гиперпараметра.
Выводы по Заданию 2: Ни одна из теоретических моделей, известных нам, не описывает структуры нашего графа достаточно точно.
В этом задании нам предлагается посчитать различные меры центральности вершин и сделать содержательные выводы на основе полученных расчетов.
Замечание: Для анализа используем обратные веса, т.к. их использование делает осмысленным расчет центральности по близости и кратчайшему пути, но не искажает остальных мер центральности вершин.
Для начала посчитаем центральность по степени:
deg <- degree(net)
# plot(net, vertex.label=NA, layout=layout_with_fr, vertex.size=deg)
Теперь построим гистограмму распределения степеней вершин:
# гистограмма степеней вершин
hist(deg, main="Histogram of node degree", breaks=50, prob=TRUE)
# распределение степеней вершин
mean(deg) # Средняя степень вершины по графу
[1] 5.655857
# ~ 5.66
Mode(deg) # Модальная степень вершин графа
[1] 1
attr(,"freq")
[1] 599
# 1 -> 599 раз
median(deg) # Медиана степеней вершин
[1] 2
# ~ 2
top_10_degree <- tail(sort(deg), 10) # топ-10 по степеням
top_10_degree # в топе - 7 университетов, 3 гос. учреждения
Canterbury District Health Board AgResearch University of Waikato
121 126 136
Auckland District Health Board AUT University Victoria University of Wellington
143 145 204
University of Canterbury Massey University University of Otago
241 310 406
University of Auckland
551
Вывод: Центральность вершин по степени не зависит от
весов. Подавляющее большинство вершин имеет степень, меньшую
50. Так, медиана кол-ва соседей равняется 2. В среднем у
вершины такого графа степень равна 5.66. Наиболее часто встречающаяся
степень вершины — 1, которая встречается 599 раз. Если мы посмотрим 10
вершин с наибольшим числом соседей, то окажется, что 7 из них являются
университетами, в числе которых крупнейший — университет Окленда
(степень вершины максимальная — 551). Оставшиеся 3 организации являются
государственными: 2 департамента здравоохранения и 1 исследовательская
компания. Степень у всех вершин, входящих в топ-10, > 100.
Далее рассчитаем центральность по близости двумя
способами: В первом варианте все веса ребер графа являются
1, в другом используются обратные веса ребер.
clos1 <- closeness(net, weights=NA) # Первый способ - без учета весов
# clos1
plot(net, vertex.label=NA, layout=layout_with_fr, vertex.size=clos1 * 50)
top_10_clos1 <- tail(sort(clos1), 10)
# Вектор-топ 42 также будет единичным, а вот вектор-топ 43 уже нет
top_10_clos1
Reveal Infrastructure Ltd. Napier City Council RATO Natural Health Clinic
1 1 1
Stahlton Engineered Concrete Katipo Communications Workwise Employment Ltd.
1 1 1
Northpower New Zealand Stock Exchange Real Time Genomics
1 1 1
Rotary Club of Riccarton
1
bottom_10_clos1 <- head(sort(clos1), 10)
bottom_10_clos1
Tauranga City Council
0.0001571092
Project Litefoot Trust
0.0001703868
Dumpark Data Science
0.0001706193
Dairy Production Systems Ltd.
0.0001720282
SciCon Scientific Consultants
0.0001725923
New Zealand Employment Service and Work and Income New Zealand
0.0001773050
IWM Consultancy
0.0001779676
PLANTwise Services Ltd
0.0001779676
Care of Orca Research Trust
0.0001794044
Southern Veterinary Centre
0.0001818843
clos1["University of Auckland"]
University of Auckland
0.0004175365
clos2 <- closeness(net)
# clos2
plot(net, vertex.label=NA, layout=layout_with_fr, vertex.size=clos2 * 50)
top_10_clos2 <- tail(sort(clos2), 10)
top_10_clos2
Workwise Employment Ltd. Northpower New Zealand Stock Exchange
1.0 1.0 1.0
Real Time Genomics Rotary Club of Riccarton Bio Soil and Crop Ltd
1.0 1.0 1.2
PollenPlus Ltd. GroPlus Ltd. Blues Rugby Club
1.2 1.5 2.0
Highlanders Rugby
2.0
bottom_10_clos2 <- head(sort(clos2), 10)
bottom_10_clos2
Dumpark Data Science
0.0002501171
Dairy Production Systems Ltd.
0.0002502698
New Zealand Employment Service and Work and Income New Zealand
0.0002504864
Project Litefoot Trust
0.0002504870
Tauranga City Council
0.0002955979
Port Otago Ltd
0.0003060701
Travel Clinic
0.0003065386
Dobbie Engineers Ltd.
0.0003069975
NZ Strong Construction
0.0003069975
Pulp Mill Technical Group
0.0003069975
clos2["University of Auckland"]
University of Auckland
0.0009338331
Вывод: В обоих случаях максимальная центральность
вершин будет составлять 1. В нашем графе существует 21 пара
вершин, которая связана лишь между собой и более ни с кем (из рассчета
для графа с единичными весами ребер). Из них лишь 5 организаций работали
совместно более 1 раза (из рассчета для графа с обратными весами).
Топ-10 центральных вершин по близости и топ-10 по степени уже
существенно различаются. Так, центральность по близости университета
Окленда составит 4.1753653^{-4} в первом случае и 4.1753653^{-4}. Т.е.
round(1 / clos1["University of Auckland"]) или
round(1 / clos2["University of Auckland"]) — сумма длин
всех путей от университета Окленда до любой другой вершины в зависимости
от метода рассчета.
betw1 <- betweenness(net, weights=NA)
head(betw1)
45 South Management Ltd. Abacus ALS
0.0000 0.0000
AbacusBio Ltd. ABB Power Quality
228.0179 0.0000
ABI Rehabilitation New Zealand Ltd. AC Research Associates
0.0000 0.0000
top_10_betw1 <- tail(sort(betw1), 10)
top_10_betw1
Auckland District Health Board Plant and Food Research AgResearch
44138.73 55148.49 61204.78
University of Waikato AUT University Victoria University of Wellington
62949.44 64866.37 84858.46
University of Canterbury Massey University University of Otago
130456.08 194334.96 231358.36
University of Auckland
385421.47
bottom_10_betw1 <- head(sort(betw1), 10)
bottom_10_betw1
45 South Management Ltd. Abacus ALS
0 0
ABB Power Quality ABI Rehabilitation New Zealand Ltd.
0 0
AC Research Associates Academic Colleges Group New Zealand
0 0
Action on Smoking and Health New Zealand Active Health
0 0
Ade Consultants Adidas Sports Medicine
0 0
betw1["University of Auckland"]
University of Auckland
385421.5
betw2 <- betweenness(net)
head(betw2)
45 South Management Ltd. Abacus ALS
0 0
AbacusBio Ltd. ABB Power Quality
1 0
ABI Rehabilitation New Zealand Ltd. AC Research Associates
0 0
top_10_betw2 <- tail(sort(betw2), 10)
top_10_betw2
Auckland District Health Board Plant and Food Research AUT University
48988.50 57655.33 74583.33
University of Waikato AgResearch Victoria University of Wellington
81434.50 96575.50 139495.17
University of Canterbury Massey University University of Otago
198521.67 323919.17 340995.17
University of Auckland
815924.25
bottom_10_betw2 <- head(sort(betw2), 10)
bottom_10_betw2
45 South Management Ltd. Abacus ALS
0 0
ABB Power Quality ABI Rehabilitation New Zealand Ltd.
0 0
AC Research Associates Academic Colleges Group New Zealand
0 0
Accident Compensation Corporation Action on Smoking and Health New Zealand
0 0
Active Health Acurity Health Group ltd.
0 0
betw2["University of Auckland"]
University of Auckland
815924.2
Вывод: При разных методах минимальные значения
центральности у вершин будут равны 0. Это вершины,
расположенные на периферии графа, или те, которые вовсе имеют всего
одного соседа (например, та 21 пара вершин). Максимальные центральности
вершин превышают 40 000 в топ-10. В обоих случаях
университет Окленда на первом месте (3.8542147^{5} и 8.1592425^{5}).
Значительное число кратчайших путей между любыми двумя вершинами
проходит именно через данный университет.
eig1 <- eigen_centrality(net, weights = NA)
# head(eig1)
top_10_eig1 <- tail(sort(eig1$vector), 10)
top_10_eig1
AgResearch Auckland District Health Board Canterbury District Health Board
0.3760716 0.4334652 0.4357203
University of Waikato AUT University Victoria University of Wellington
0.4373508 0.4700463 0.5846012
University of Canterbury Massey University University of Otago
0.5943213 0.6779365 0.8521192
University of Auckland
1.0000000
bottom_10_eig1 <- head(sort(eig1$vector), 10)
bottom_10_eig1
Highlanders Rugby Growth Solutionz
3.407685e-18 3.423996e-18
FEM Ltd. New Zealand Rock Lobster Industry Council (NZ RLIC)
3.524028e-18 3.563296e-18
GHD Limited Napier City Council
3.577456e-18 3.578251e-18
Rotary Club of Riccarton Milner Consulting Ltd.
3.585397e-18 3.585975e-18
NetValue Crown Fibre Holdings Limited
3.598880e-18 3.629860e-18
eig2 <- eigen_centrality(net)
# head(eig2)
top_10_eig2 <- tail(sort(eig2$vector),10)
top_10_eig2
NIWA University of Waikato AUT University
0.1576988 0.1897041 0.1992335
Canterbury District Health Board Auckland District Health Board Victoria University of Wellington
0.2127479 0.2346215 0.2896466
University of Canterbury Massey University University of Otago
0.3085040 0.4043171 0.6707261
University of Auckland
1.0000000
bottom_10_eig2 <- head(sort(eig2$vector),10)
bottom_10_eig2
Alzheimers Northland AP EnerTec Ltd.
0 0
Aviation Industry Association of New Zealand Inc Bio Soil and Crop Ltd
0 0
Blues Rugby Club Bruce Shallard and Associates
0 0
ChangeMakers Refugee Forum Chatham Rock Phosphate Ltd.
0 0
Chilton Saint James School Crown Fibre Holdings Limited
0 0
Вывод: Минимальные значения центральностей по
собственному значению также составляют 0. Наибольшую же
имеет университет Окленда. Его высокий показатель говорит, что сам
университет связан со многими вершинами, которые в свою очередь имеют
также высокое значение центральности.
# объединим все центральности в один датафрейм и посмотрим корреляции в двух вариантах
df1 <- data.frame(deg, clos1, betw1, eig1$vector)
cor(df1)
deg clos1 betw1 eig1.vector
deg 1.00000000 -0.03452649 0.96148582 0.8926645
clos1 -0.03452649 1.00000000 -0.01592647 -0.1091064
betw1 0.96148582 -0.01592647 1.00000000 0.7596773
eig1.vector 0.89266448 -0.10910642 0.75967732 1.0000000
df2 <- data.frame(deg, clos2, betw2, eig2$vector)
cor(df2)
deg clos2 betw2 eig2.vector
deg 1.00000000 -0.03349182 0.90929438 0.8418800
clos2 -0.03349182 1.00000000 -0.01198628 -0.1332681
betw2 0.90929438 -0.01198628 1.00000000 0.7737961
eig2.vector 0.84188004 -0.13326806 0.77379607 1.0000000
Вывод: Составив корреляционные матрицы заметим, что существует высокая корреляция между следующими парами центральностей вершин: по степени - по кратчайшему расстоянию, по степени - по собств. значению и по собств. значению - по кратчайшему расстоянию
Сравним, какие организации попадали в топы по центральностям (исключая центральность по близости).
corr1 <- data.frame(
names(top_10_degree),
names(top_10_betw1),
names(top_10_eig1)
)
corr1
corr2 <- data.frame(
names(top_10_degree),
names(top_10_betw2),
names(top_10_eig2)
)
corr2
Выводы по Заданию 3: Сравнив топы по максимальным центральностям (за исключением топа центральностей по близости из-за малого количества соседей [порядка 1]), заметим, что по трем центральностям первую пятерку занимают одни и те же организации: университет Окленда, университет Отаго, университет Мэсси, университет Кантербери и Викторианский университет Веллингтона.
В Задании 4 нам предлагается рассчитать ассортативность и гомофилию для графа. Начнем исследование данных с ассортативности:
# Возникли проблемы с вычислением ассортативности через
# assortativity_nominal(net, types=V(net)$type, directed=F)
# Предупреждение: в результате преобразования созданы NA
# [1] NaN
assortativity_nominal(net, types=V(net)$type, directed=F)
Предупреждение: в результате преобразования созданы NA
[1] NaN
# Поэтому было принято решение создать новый столбец с закодированными строками
V(net)$num_type <- ifelse( V(net)$type=="Business Enterprise", 1,
ifelse( V(net)$type=="Private not for profit", 2,
ifelse( V(net)$type=="Government", 3,
ifelse( V(net)$type=="Higher Education", 4, 0)
)
)
)
assort_nom <- assortativity_nominal(net, types=V(net)$num_type, directed=F)
assort_nom # -0.05876298
[1] -0.05876298
# значение ассортативности не зависит от конкретных чисел
V(net)$num_type1 <- ifelse( V(net)$type=="Business Enterprise", 2,
ifelse( V(net)$type=="Private not for profit", 7,
ifelse( V(net)$type=="Government", 5,
ifelse( V(net)$type=="Higher Education", 10, 0)
)
)
)
assortativity_nominal(net, types=V(net)$num_type1, directed=F) # -0.05876298
[1] -0.05876298
# считаем ассортативность по степени вершины
assort_deg <- assortativity_degree(net, directed=F) # -0.3361749
assort_deg
[1] -0.3361749
Обобщение результатов: Ассортативность в графе по типу вершин примерно равна -0.059. Это означает, что мы не можем утверждать наличие предпочтительного присоединения узлов сети к своему типу (среди узлов в общем, т.е. среди узлов всех типов).
Ассортативность по степени вершины составляет примерно -0.336, что говорит о предпочтительном присоединении узлов к вершинам с относительно небольшим количеством соседей.
Для выявления более точного эффекта предпочтительного присоединения именно среди университетов (что интересует нас в рамках исследовательского вопроса) рассчитаем показатели гомофилии для университетов и других типов институтов.
Нами было принято решение рассчитать dyadicity (ее посчитаем
только для одного случая, т.к. она не будет меняться при различных
разбиениях) и heterophilicity учреждений высшего образования
для следующих разделений на подгруппы: - Higher Education
vs все остальные - Higher Education vs
Business Enterprise - Higher Education vs
Government.
Замечание: Влиянием Private not for profit можно
пренебречь.
Higher Education vs все остальныеp <- edge_density(net)
p # 0.003745601
[1] 0.003745601
n_c <- sum( V(net)$type == 'Higher Education' )
n_c # 54
[1] 54
exp_diad_uni <- n_c * (n_c - 1) / 2 * p
exp_diad_uni # 5.359955
[1] 5.359955
unis <- V(net)[which( V(net)$type == 'Higher Education' )]
unis <- names(unis)
# unis
m_c <- 0
for (row in 1:nrow(edges)) {
from_obj <- edges[row, 'from']
to_obj <- edges[row, 'to']
if ( (from_obj %in% unis) && (to_obj %in% unis) ) {
m_c <- m_c + 1
}
}
m_c # 172
[1] 172
diadicity_uni <- m_c / exp_diad_uni
diadicity_uni # 32.08982
[1] 32.08982
# 2-й способ посчитать dyadicity
uni_net <- induced_subgraph(
net,
v=V(net)[which( V(net)$type == 'Higher Education' )]
)
# uni_net
# V(uni_net)$type # проверяем, что все вершины в новом графе вузы
edge_density(net) # 0.003745601
[1] 0.003745601
edge_density(uni_net) # 0.1201957
[1] 0.1201957
edge_density(uni_net) / edge_density(net) # 32.08982
[1] 32.08982
n_c # 54
[1] 54
p # 0.003745601
[1] 0.003745601
n_nc <- sum( V(net)$type != 'Higher Education' )
n_nc # 1457
[1] 1457
exp_hetero_uni <- n_nc * n_c * p
exp_hetero_uni # 294.6964
[1] 294.6964
m_mixed <- 0
for (row in 1:nrow(edges)) {
from_obj <- edges[row, 'from']
to_obj <- edges[row, 'to']
if ( (from_obj %in% unis) && !(to_obj %in% unis) ) {
m_mixed <- m_mixed + 1
}
if ( !(from_obj %in% unis) && (to_obj %in% unis) ) {
m_mixed <- m_mixed + 1
}
}
m_mixed # 2051
[1] 2051
heterophilicity_uni <- m_mixed / exp_hetero_uni
heterophilicity_uni # 6.959706
[1] 6.959706
Higher Education vs
Business Enterprise# uni_biz_net <- induced_subgraph(
# net,
# v=V(net)[which(
# (V(net)$type == 'Higher Education') || (V(net)$type == 'Business Enterprise')
# )]
# )
uni_biz_net <- induced_subgraph(
net,
v=V(net)[which( V(net)$type %in% c('Higher Education', 'Business Enterprise') )]
)
# uni_biz_net
# V(uni_biz_net)$type
uni_wout_biz_net <- induced_subgraph(
uni_biz_net,
v=V(uni_biz_net)[which( V(uni_biz_net)$type == 'Higher Education' )]
)
# uni_wout_biz_net
# V(uni_wout_biz_net)$type
p_uni_vs_biz <- edge_density(uni_biz_net)
p_uni_vs_biz # 0.003036165
[1] 0.003036165
edge_density(uni_wout_biz_net) # 0.1201957
[1] 0.1201957
dyadicity_uni_vs_biz <- edge_density(uni_wout_biz_net) / p_uni_vs_biz
dyadicity_uni_vs_biz # 39.58799
[1] 39.58799
n_c # 54
[1] 54
p_uni_vs_biz # 0.003036165
[1] 0.003036165
n_nc_biz <- sum( V(net)$type == 'Business Enterprise' )
n_nc_biz # 1010
[1] 1010
exp_hetero_uni_vs_biz <- n_nc_biz * n_c * p_uni_vs_biz
exp_hetero_uni_vs_biz # 165.5925
[1] 165.5925
biz <- V(net)[which( V(net)$type == 'Business Enterprise' )]
biz <- names(biz)
# biz
m_mixed_vs_biz <- 0
for (row in 1:nrow(edges)) {
from_obj <- edges[row, 'from']
to_obj <- edges[row, 'to']
if ( (from_obj %in% unis) && (to_obj %in% biz) ) {
m_mixed_vs_biz <- m_mixed_vs_biz + 1
}
if ( (from_obj %in% biz) && (to_obj %in% unis) ) {
m_mixed_vs_biz <- m_mixed_vs_biz + 1
}
}
m_mixed_vs_biz # 1139
[1] 1139
heterophilicity_uni_vs_biz <- m_mixed_vs_biz / exp_hetero_uni_vs_biz
heterophilicity_uni_vs_biz # 6.878333
[1] 6.878333
Higher Education vs Governmentuni_gov_net <- induced_subgraph(
net,
v=V(net)[which( V(net)$type %in% c('Higher Education', 'Government') )]
)
# V(uni_gov_net)$type
uni_wout_gov_net <- induced_subgraph(
uni_gov_net,
v=V(uni_gov_net)[which( V(uni_gov_net)$type == 'Higher Education' )]
)
# V(uni_wout_gov_net)$type
p_uni_vs_gov <- edge_density(uni_gov_net)
p_uni_vs_gov # 0.04231381
[1] 0.04231381
edge_density(uni_wout_gov_net) # 0.1201957
[1] 0.1201957
dyadicity_uni_vs_gov <- edge_density(uni_wout_gov_net) / p_uni_vs_gov
dyadicity_uni_vs_gov # 2.840578
[1] 2.840578
n_c # 54
[1] 54
p_uni_vs_gov # 0.04231381
[1] 0.04231381
n_nc_gov <- sum( V(net)$type == 'Business Enterprise' )
n_nc_gov # 1010
[1] 1010
exp_hetero_uni_vs_gov <- n_nc_gov * n_c * p_uni_vs_gov
exp_hetero_uni_vs_gov # 2307.795
[1] 2307.795
gov <- V(net)[which( V(net)$type == 'Business Enterprise' )]
gov <- names(gov)
# gov
m_mixed_vs_gov <- 0
for (row in 1:nrow(edges)) {
from_obj <- edges[row, 'from']
to_obj <- edges[row, 'to']
if ( (from_obj %in% unis) && (to_obj %in% gov) ) {
m_mixed_vs_gov <- m_mixed_vs_gov + 1
}
if ( (from_obj %in% gov) && (to_obj %in% unis) ) {
m_mixed_vs_gov <- m_mixed_vs_gov + 1
}
}
m_mixed_vs_gov # 0.4935447
[1] 1139
heterophilicity_uni_vs_gov <- m_mixed_vs_gov / exp_hetero_uni_vs_gov
heterophilicity_uni_vs_gov # 0.4935447
[1] 0.4935447
Соберем полученные значения dyadicity и
heterophilicity вместе: - Higher Education vs все
остальные: - dyadicity: 32.0898241 - heterophilicity: 6.9597056 -
Higher Education vs Business Enterprise: -
dyadicity: 39.5879866 - heterophilicity: 6.8783329 -
Higher Education vs Government: - dyadicity:
2.840578 - heterophilicity: 0.4935447.
Мы видим, что dyadicity университетов при разделении на университеты и всех остальных (а также при разделении на университеты и бизнес) очень высокая. При этом heterophilicity также доволльно сильно превышает 1.
В чем может быть объяснение этой особенности данных? Как утверждается
в статье “Functional
dyadicity and heterophilicity of gene-gene interactions in statistical
epistasis networks”, описанный выше факт может быть связан с тем,
что вершины со значением метки 1 имеют большую
центральность \(\Rightarrow\) хорошо
связаны друг с другом и с вершинами другого типа.
В Задании 6 мы будем проверять гипотезу касательно высоких
степеней вершин типа Higher Education среди крупных
университетов.
Выводы по Заданию 4: Мы получили, что несмотря на плотные связи учреждений высшего образования между собой, по крайней мере некоторые из них имеют также прочные связи с другими типами организаций (в частности с бизнесом).
Проведем процедуру кластеризации вершин при помощи Edge-betweenness и Fast-greedy методов и сравним полученные результаты с существующей классификацией.
Т.к. считать кратчайший путь для исходных весов ребер является не очень осмысленным, используем обратные для исходных веса ребер.
Замечание: Для Fast-greedy не важно, какие данные использовать: исходные или инвертированные, т.к. этот метод не учитывает веса ребер (только их количество).
# проводим кластеризацию при помощи метода Edge-betweenness
ceb <- cluster_edge_betweenness(net, weights=E(net)$weight)
Предупреждение: At core/community/edge_betweenness.c:493 : Membership vector will be selected based on the highest modularity score.
# смотрим на полученные в результате кластеризации объекты
mbship_ceb <- membership(ceb)
# mbship_ceb
s_ceb <- sizes(ceb)
s_ceb
Community sizes
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26
45 213 233 154 229 190 75 46 164 48 32 2 2 2 3 15 2 2 8 2 2 4 2 2 2 2
27 28 29 30 31 32 33 34 35 36 37 38
7 2 2 2 2 2 2 2 2 2 3 2
ms_ceb <- modularity(ceb)
ms_ceb # 0.527262
[1] 0.527262
# ceb
# > s_ceb <- sizes(ceb)
# > s_ceb
# Community sizes
# 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
# 45 213 233 154 229 190 75 46 164 48 32 2 2 2 3 15 2 2 8 2
# 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38
# 2 4 2 2 2 2 7 2 2 2 2 2 2 2 2 2 3 2
# > ms_ceb <- modularity(ceb)
# > ms_ceb
# [1] 0.527262
Визуализация кластеризации методом Edge-betweenness:
plot(
ceb, net,
vertex.label=NA,
vertex.color=V(net)$color,
layout=layout_nicely
)
plot(
ceb, net,
vertex.label=NA,
vertex.color=V(net)$color,
layout=layout_with_lgl
)
Сравним результаты кластеризации Edge-betweenness и имеющуюся классификацию вершин:
compare(ceb, V(net)$type, method='nmi') #
[1] 0.03399001
# здесь и далее мы используем нормализованное значение функции compare() для облегчения интерпретации возвращаемого результата
Таким образом, мы видим, что имеющаяся классификация вершин совсем не похожа на то, что мы получили при кластеризации. По всей видимости, “разделение по бутылочному горлышку” является не очень хорошим способом для кластеризации вершин в нашем датасете.
# проводим кластеризацию при помощи метода Fast-greedy
cfg <- cluster_fast_greedy(net, weights=E(net)$weight)
# смотрим на полученные в результате кластеризации объекты
mbship_cfg <- membership(cfg)
# mbship_cfg
s_cfg <- sizes(cfg)
s_cfg
Community sizes
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26
219 144 212 132 141 94 189 120 84 82 28 3 2 15 3 3 2 2 2 2 2 2 2 2 2 2
27 28 29 30 31 32 33 34 35 36
2 2 2 2 2 2 2 2 2 2
ms_cfg <- modularity(cfg)
ms_cfg # 0.5310055
[1] 0.5310055
# cfg
# > s_cfg <- sizes(cfg)
# > s_cfg
# Community sizes
# 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26
# 219 144 212 132 141 94 189 120 84 82 28 3 2 15 3 3 2 2 2 2 2 2 2 2 2 2
# 27 28 29 30 31 32 33 34 35 36
# 2 2 2 2 2 2 2 2 2 2
# > ms_cfg <- modularity(cfg)
# > ms_cfg
# [1] 0.5310055
Визуализация кластеризации методом Fast-greedy:
plot(
cfg, net,
vertex.label=NA,
vertex.color=V(net)$color,
layout=layout_nicely
)
plot(
cfg, net,
vertex.label=NA,
vertex.color=V(net)$color,
layout=layout_with_lgl
)
compare(cfg, V(net)$type, method='nmi') # 0.03386511
[1] 0.03386511
Таким образом, мы видим, что так же, как и предыдущий способ кластеризации, Fast-greedy не очень хорошо соответствует реальной классификации рассматриваемых институтов.
Вывод по Заданию 5: Известные нам способы кластеризации (Edge-betweenness и Fast-greedy) не дали результата, похожего на имеющуюся в датасете классификацию.
Изначально был поставлен исследовательский вопрос, насколько активно при написании научных работ университеты Новой Зеландии взаимодействуют с другими типами организаций, нужно ли им как-то в этом способстовать извне (например, создавать общие лаботарии с индустриальными корпорациями и государственными органами). Помимо расчета значений гомофилии и различных видов центральностей, для ответа на заданный вопрос нам может помочь модель предсказания класса вершин Relational Neighbors Classifier (RNC).
Разделим все вершины на 2 класса: “1” - университет, “0” - иначе.
Отберем из каждого класса по 10 вершин и положим вероятности их
принадлежности первому классу соответственно 1 и
0. Для всех остальных вершин (unknown) положим
эту вероятность равной 0.5. При каждой итерации будем
обновлять вероятности каждой вершины исходя из вероятностей её соседей.
В итоге, после 10 таких итераций получим для каждой вершины вероятность
её принадлежности классу “1”, и отберем те вершины, чья вероятность
больше порога 0.5 в качестве принадлежащих классу “1”.
Для избежания случайности повторим этот алгоритм 100 раз и построим гистограмму распределения числа предсказанных университетов, а также гистограмму распределения средних вероятностей принадлежности группе университетов.
res <- rep(NA, 20)
pr <- rep(NA, 20)
for (j in 1:20){
network <- net
V(network)[which( V(network)$type == 'Higher Education' )]$type = 1
V(network)[which(V(network)$type != 1)]$type = 0
# вероятность вершины перейти в класс churn = 1 равна доле соседей churn = 1 от общего числа соседей
# используем матричные вычисления
random_uni <- sample( V(network)[ which(V(network)$type == 1) ], 10 )
random_not_uni <- sample( V(network)[ which(V(network)$type == 0) ], 10 )
for (i in 1:1511){
if (!(V(network)[i] %in% random_uni | (V(network)[i] %in% random_not_uni))){
V(network)$type[i] = 0.5
}
}
UniProb <- as.numeric(V(network)$type)
threshold <- 0.5
A <- as_adjacency_matrix(network)
A <- as.matrix(A)
Neighbors <- rowSums(A) # общее число соседей вершины
UniProb_iter <- UniProb
for(i in 1:10){
UniProb_iter <- as.vector((A %*% UniProb_iter) / Neighbors)
}
res[j] <- length(which(UniProb_iter > threshold))
pr[j] <- mean(UniProb_iter)
}
hist(res, main='Число предсказанных университетов')
hist(pr, main='Средняя вероятность принадлежности классу 1')
Как можно увидеть, почти во всех случаях абсолютное число вершин предсказанj в качестве университетов. Получаем, что какие бы 10 университетов мы ни взяли, их класс дойдет до большинства вершин быстрее, чем противоположный класс 10 случайных вершин. Это распространяется и на случай, когда были выбраны наименее центральные университеты. Если бы университеты были плохо связаны с организациями другого типа, то их класс намного слабее распространялся бы на вершины другого типа, а это не так.
Значит, по крайней мере какие-то университеты (а именно наиболее центральные) довольно активно кооперируются с другими организациями (как, например, это делает Оклендский университет). Поэтому для таких университетов нет необходимости строить совместные с бизнесом и государством научные центры.
Вспомним про высокий показатель dyadicity у ВУЗов: выходит, что он получается таким большим из-за высокой гомофилии более мелких учебных заведений. Тогда можно предположить, что именно они плохо взаимодействуют с бизнесом и гос. органами.
Как тогда можно им помочь? При моделировании их класс распространялся по сети через их соседей — крупные университеты. В таком случае, более разумно строить научные центры именно соединяющие мелкие университеты с крупными. Данным образом, казалось бы, усиливая гомофилию, мы, наоборот, будем способстовать развитию научных отношений университетов с индустриальными партнерами.
Здесь мы предполагаем эталонность крупных университетов (например, Оклендского) с точки зрения эффективности взаимодействия, а также качества и количества научных работ, что продемонстрировано, например, в этой статье.